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Abstract 

We calculate the relativistic correction to J/ip production in the color-singlet process e + e~ —> 
J/i/j + gg at l?-factories. We employ the nonrelativistic QCD factorization approach, where the 
short-distance coefficients are calculated perturbatively and the long-distance matrix elements are 
extracted from the decays of J/ip into e + e~ and light hadrons. We find that the 0(v 2 ) relativistic 
correction can enhance the cross section by a factor of 20-30%, comparable to the enhancement 
due to the 0(a s ) radiative correction obtained earlier. Combining the relativistic correction with 
the QCD radiative correction, we find that the color-singlet contribution to e + e~ —> J/ip + gg can 
saturate the latest observed cross section a(e + e~ — > J/ifi+X aon - C c) = 0.43±0.09±0.09 pb by Belle, 
thus leaving little room to the color-octet contributions. This gives a very stringent constraint on 
the color-octet contribution, and may imply that the values of color-octet matrix elements are 
much smaller than expected earlier by using the naive velocity scaling rules or extracted from 
fitting experimental data with the leading-order calculations. 
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I. INTRODUCTION 



The nonrelativistic QCD (NRQCD) effective field theory, introduced in f 995 by Bodwin, 
Braaten, and Lepage [lj], has been widely accepted as a rigorous approach to study the 
production and decay of heavy quarkonium, the bound state of heavy quark Q and antiquark 
Q pair (see Ref.j^] for a review). In the framework of NRQCD the production of heavy 
quarkonium is factorized into two parts, the short- distance part and the long-distance part. 
In the short- distance part, the QQ pair is created in certain J PC and color states, which 
can be calculated perturbatively through the expansion of QCD coupling constant a s . The 
long-distance part describes the evolution of the QQ pair into the physical hadron states 
through the emission of soft gluons with the corresponding universal nonperturbative matrix 
elements, which are weighted by powers of the relative velocity v Q of heavy quarks in the 
meson rest frame. One intrinsic character of NRQCD is the inclusion of the effect of QQ 
pair in a color-octet state, i.e. the color-octet mechanism. Since in the e + e~ collision, the 
structure of the parton involved is simpler and the signals can be prominent, it is a good 
place to study the heavy quarkonium production and test the color-octet mechanism. 

In recent years, the cross sections of inclusive J/^ production in e + e~ annihilation at 
yfs = 10.6GeV have been reported by BaBarQ], BelleQ, and CLEoQ collaborations with 
the results 

a(e + e- J/ip + X) = 2.5 ± 0.21 ± 0.21pb, (la) 
(j(e + e" -> J/ip + X) = 1.47 ± 0.10 ± 0.13pb, (lb) 
a(e + e- J/ip + X) = 1.9 ± 0.2pb, (lc) 



respectively. If the J/^'s momentum p*^ is restricted to p*^ > 2GeV/c, the results of BaBar 
and Belle |4| become 

a(e + e- M + X) = 1.87 ± 0.10 ± 0.15pb, (2a) 

P ;>2GeV/c 

ff(eV^JM + I) = 1.05 ±0.04 ±0.09pb. (2b) 

P ;>2GeV/c 

It can be seen that these experimental measurements are not entirely consistent with each 
other. The two I?-factories also measured the J/ip momentum distributions and observed 



similar shapes. Despite of the disagreement among different experimental measurements, 
theoretically the inclusive J / ip production has been extensively investigated within the color- 



singlet model[6|-|9( and the color-octet mechanism in NRQCD [9b11|]. In NRQCD, one impor 



tant contribution to the inclusive J/ip production in e + e annihilation at ^/s = 10.6 GeV 
comes from the color-octet process e+e~ J/ip^S^, 3 Pj 8) }+9 BQ- But its predictions of 
the J/ip enhancement near the kinematics end point region is not observed. After applying 
the soft-collinear effective theory (SCET) 12j, the shape of J/ip momentum distribution can 
be softened. However this depends phenomenologically on a nonperturbative shape function. 

Experimentally, the e + e~ — > J/ip + X process can be divided into two parts: the e + e~ — > 
J/ip + cc + X part and the e + e~ — >■ J/ip + X non _ cE part. The Belle collaboration finds the 
ratio fl3l| 

_ a(e + e- -»■ J/ip + cc + X) 
cE ~ a{e+e- =► J/ijj + X) 



0.591^ ±0.12. 



(3) 



which corresponds to a(e + e — > J/vp + cc+X) = 0.87_ ; 19 ±0.17pb. The latest measurement 
of J/ip production in association with a cc pair carried out by Belle givesjfjj 

cr(e + e- -> J/ip + cc + X) = 0.74 ± 0.08±g;g|pb. (4) 



The experimental results are more than 5 times larger than leading-order (LO) NRQCD 

nnn 

predictionsp], |8|, H5(- This large discrepancy could be resolved by including the NLO QCD 

Belle also analyzed the 



corrections and the feed-down of higher excited states[16l. 
e + e~ — > J/ip + X non „ cE process and obtained [l^ 



1 E3. 



cT(e + e -> J/ip + X n 



0.43 ±0.09 ±0.09 pb. 



(5) 



At LO in a s the non — cc process via color-singlet channels on 



y includes e + e — >■ J/ib + 



18 



19] 



gg. The theoretical prediction at LO in a s and v 2 is about 0.2pb[9|, and recent works 
find that the K-f actor of the NLO correction is about 1.2-1.3. In charmonium system 
the relative velocity v c of a cc pair in the J/ip rest frame is not small. The value of v 2 
is about 0.3, which is close to the size of a s (m c ). So the relativistic effect may also be 



important. Our previous work 20| shows that in the e + e~ — > J/ip + r] c exclusive process, the 
relativistic correction is indeed important in resolving the more than one order of magnitude 
discrepancies between experimental data 2l|, |22J and LO NRQCD predictions 231] . We find 
when combining the relativistic correction together with the NLO QCD corrections 24J the 
conflict between experimental measurement and theoretical prediction is almost resolved (see 



also Ref. 25[ for a similar result). However, in the e + e — > J/ip + cc process, we find the 



relativistic correction is very small and may be ignored. And the recent work shows 26] that 
in the pp — > J/ip+X process at the Tevatron the relativistic corrections can also be neglected. 
Then it is necessary to investigate the relativistic correction to the e + e~ — > J/ip + gg process 
and clarify how large the relativistic correction is and whether it is positive or negative. In 
this work we will deal with this problem within the framework of NRQCD factorization 
approach, and in particular, we will examine the effect of the relativistic correction on the 
constraint on color-octet matrix elements. The rest of this paper is organized as follows. In 
Sec. II, we will introduce the NRQCD factorization formula and describe how we calculate 
the short-distance coefficients and determine the long-distance matrix elements. We will 
present our calculations and show our result in Sec. III. Discussions and a summary will be 
given in the last section. 

II. THE NRQCD FACTORIZATION FORMULA 

According to NRQCD effective theory, up to 0{v 2 ) the inclusive J/ip production rate 
can be expressed as 

a(e+e - ^ J/^ + X) = ^C^(o|0f( 8 Si)|O) + 9l^A(0\vf( 3 S 1 )\0) +0(v 4 a). (6) 
where 

Of ( 3 S0 = x Vty( a t a^V x (7) 

and 

-Pti'Si) = \[x ] <J l ^{a\a^a\- l -^)\ + *V(~ ft) V(4^V X ] (8) 

are the four-fermion operators with dimensions six and eight, respectively. F\( 3 Si) and 
G±( 3 Si) are the corresponding short- distance coefficients. The short-distance coefficients 
can be obtained perturbatively through the matching condition 

<QQ) = E ^H(0|O^(A)|0) (9) 

pcrtQCD ^ m d ™- AX n pert NRQCD 

n L 

The long-distance matrix elements characterized by the velocity v c can be estimated by 
lattice calculations or phenomenological models, or determined by fitting experimental data. 



In this work, the covariant spinor projection method 27] is adopted to evaluate the left 



hand side of Eq. [9]. In this method, the Dirac spinor product v(P/2 — q)u(P/2 + q) is 



projected onto a certain ( 2S+1 ^Lj state in a Lorentz covariant form (see, e.g., 28|, |29|), which 
makes the short-distance coefficients evaluated directly. In the J/ijj (with S — 1) case, the 
expression of the spinor production projection in the meson rest frame up to all orders of v 2 



is 



29] 



w (-q> A 2 )ii(q, Ai)(-, Ai; -, A 2 |l, e) = 

Ai A2 

1/^ w aq s 1 + 70 . c*-q N , , 

+ __£.)«. e _^ (1 + _S. h , (10) 

And in an arbitrary frame, it becomes 

Y w (<?, A 2 )w(g, Ai)(i, Ai; ^, A 2 |l, e) = 

Ai A2 

-/» - - m)i^— —{-j> + i + m). (11) 



2V2(E + my2" * >r 2E V 2 J 

Here the normalization of the Dirac spinors is uu = —vv = 2m c , and the relations between 
momenta of quark and antiquark in an arbitrary frame and in the meson rest frame are 
given by [301] 

l -P + q = L{ l -P r + c L ) ) (12a) 

\p-q = L{ l -P r -c L ), (12b) 

where Pfi = (2E q , 0), E q = \J m? + q 2 , and 2q is the relative momentum between two quarks 
in the meson rest frame. is the boost tensor from the meson rest frame to an arbitrary 
frame. 

Then the cross section for e + e _ — > J/tp + gg up to next-to-leading order in v 2 can be 
expressed as 

ff(e + e - m+m) = mpmi_ j Wam + (o^wi 1^ ^ (13) 

where $3 is the three-body phase space, bar means averaging the spins over the initial 
states and summing up the spins over the final states. The short-distance part N and 
N\ defined in Eq. |17|fT8] can be calculated perturbatively. The numerical values of the 
two long-distance matrix elements (0\Of ( 3 5i)|0) and (0\Vf ( 3 5i)|0) will be estimated by 
nonperturbative methods. In the nonrelativistic limit, (0\Of ( 3 Si)|0) can be related to the 
nonrelativistic bound state wave function at the origin. And in NRQCD effective theory, 
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FIG. 1: The typical two among six Feynman diagrams in the process e + e — > J /ip + gg at leading 
order in QCD. 



(0|pf( 3 ^)|0) also can be calculated with the help of potential model[3lJ. Here we determine 
them with an alternate way by extracting them from the decays of J/ip into light hadrons 
and into e + e~. 



III. RELATIVISTIC CORRECTIONS TO e+e" -> J/ip + gg 
A. Short-Distance Coefficients 

At LO in QCD, there are six Feynman diagrams in the process e + e~ — > J/ip + gg, and the 
typical two are shown in Fig.[l]. The full QCD amplitude for this process can be expanded 
in terms of the quark relative momentum q^: 

M (e + e- (CC)a 5l (P J/V) ) +gg) = Q) 1,2 A{q^) = 



w„, „dA 



where 



1 „ r d 2 A 



^ {m+ ^_^ + .^ (14) 



q-4, =o 



A fa) = E E^' Al ' l> As|1 ' ^ )(3 ' Z ' ; 3 ^\l)A(e + e 



'2 2 

AiA 2 ij 



(15) 



where (3,2;3,j|l) = Sij/y/Nc is the color-SU(3) Clebsch-Gordon coefficient for a cc pair 

projecting onto a color-singlet state. With the help of Eq.pXj, we can express A(q^) in a 

covariant form. The factor (— ^-) 1 ^ 2 comes from the relativistic normalization of the cc state, 

E 

and E = ^m? c + q^ 2 . 



For the S wave state the odd-power terms of qj/^ vanish and q a q 13 = ^-{—g af5 + pa p P ^ ) = 
^n a/3 , where P 2 = 4E 2 , P ■ q = 0. Then at LO of v 2 we have a 

\M\ 2 = ^A(0)A*(0) + \q^A a pA\<S) + 1^4^4(0), (16) 

where A a p = 9 ®Iq p , and A* a ^ = d d *J£ p ■ According to the spinor projection method the 
short-distance part of the v® part is 

No = w ^—(A(0)A*(0)) , , (17) 



2N c m c 



and the corresponding v 2 part is 
_ = 1 t d(rfA(0)A*(0)) 



2N c m c K 9(qJ) 



1 

o ' 6 



+ ^(AfA*(0) + A(0)A^) , ). (18) 



1.7/ V " 

We introduce the dimensionless variables Zi = 2£ , j/v / s, = 2pl / ^/s,xi = cos^j and 
5 = Am c / yfs to describe the e~(k 1 ) + e + (A; 2 ) — > J/ip(p 1 ) + (7(^2) + g{Pz) process. Here y/s 
is the total energy in the center of mass frame, pi, p 2 , p 3 are the four-momenta of the 
final state J/ip and the two gluons respectively, and 6i is the angle between state i and the 
electron. The scalar products between the momenta can be expressed as 

S S S 

h-Pi = j(zi - qiXx); k 2 • Pi = -{z 1 + qiXi); p 2 ■ p 3 = -(A - Az l + 5 2 ); 

s s 
h-P2 = -^{z 2 - (q-x_ - q 1 x 1 )/2); k 2 ■ p 2 = ~(z 2 + (q~x_ - q^/2); 

S S 

h-P3 = ^(23 + (q-x- + q 1 x 1 )/2); k 2 ■ p 3 = -(z 3 - (g_x_ + q^^/2); 
Pi-P2= S -(4-5 2 - 4z 3 ); Pl .p 3 = S -(A-5 2 - 4z 2 ); h ■ k 2 = £; (19) 

where z_ = z 2 — z 3 , q_ = \q 2 — q 3 \ = a/4 — Az[ + 5 2 + z 2 _ , q\ — \q\\ — \Jz\ — 8 2 , X- = cos#_, 
and 6L is the angle between q_ and the electron. And the three-body phase space is then 
given by 

3 

2 E, 



d<$> 3 = (2vr) 4 5\h + k 2 - Pl -p 2 - p 3 ) H < '"' J ' 



i=i 

s dz\dx\dz-dw 
32(2tt) 4 yfty- K 2 )(l - xj) 



(20) 



We do not expand the three-body phase space by v 2 c and assume Mj/^ = 2m c . 
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where 

K = Z - (2 - Z ^\ (21a) 

qiq- 

w = X-+Kx x . (21b) 
The ranges of those integral variables are 

5 2 

8 <z x <l + (22a) 
- 1 < xi < 1, (22b) 



-J(zl-5*)<z^<J(zf-S2), (22c) 



(1 - K 2 )(l - xj) < w < \ {l - K 2 )(l - x 2 ). (22d) 



B. Long-Distance Matrix Elements 

The color-singlet production matrix elements can be related to the decay matrix elements 
in the vacuum saturation approximation [1] and the errors are of v\ order. Since in this 
e + e~ — > J/ip + gg process, there are two NRQCD matrix elements accurate to order v 2 , 
i.e., (0|(9f ( 3 Si)|0) and (0\Vf ( 3 «Si)|0), we can determine their values by fitting J/i/j decays 



into e + e and into light hadrons 
J/ip -> e + e" and J/ip -> LH b are 



H). The theoretical results at NLO in a* and v 2 for 



29| 



r(.W - e^) = »(1 - ^)i^MZ? _ (23 a) 

3 V 37r 3 / 

v / ; v 243 V ;; V V tt ; m 2 r 12tt 2 -108 m\ J 

(23b) 



And the central values of the experimental results are 



32] 



r(J/V> -> e+e") = 5.54keV, T(J/ip -> Li/) = 69.2keV (24) 
Solving these equations at LO of a s (QCD radiative corrections not included), we obtain 
mU^B = o. 294 GeV», f» = 0.320 x 10-GeV, (25) 



3m? 



c 



b We do not include the electromagnetic process J/-0 — > 7* — > 
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for m c =1.5GeV and a s = 0.26. Fixing a s = 0.26, we can express the matrix elements as 
functions of m c , which are 

m^lM = . 131m? , <°l^( 3 f')l°> = 0.142 x 10-V- (26) 
3 c 3m 2 c K J 

If we include the QCD NLO radiative corrections in Eq. [23], for m c =1.5Gev and a s = 0.26, 
we get 

■)Vv3c Mn\ /nlTWY3i 



mmm = . 572GeV3 , mm = . 512 x 10 - I(3eV . 

3 3m;: 
The corresponding m c dependence of the matrix elements are 

m^m = . 254m? , MM = . 228 x 10 - v . 

3 3m~ 



(27) 



(28) 



In Ref. 



33|, the authors relate (0\pf ( 3 Si)|0) to (0|Cf ( 3 Si)|0) using the NRQCD equation 



of motion . In Ref. 



3l{|, {0\Vf ( 3 Si)|0) is calculated by combining NRQCD and the potential 



model. Alternately, we estimate the production matrix elements by using the experimentally 
observed charmonium decay rates. As argued in[l|, the differences between the color-singlet 
production and decay matrix elements are of order t> 4 . So our method should be valid at 
order v 2 , and our numerical results of the production matrix elements are adequate, if the 
higher order QCD and v 2 corrections to the decay rates are small and the uncertainties of 
the experimental data are not large. If we define 0((v 2 )) = ^^Lag^ , we find the value 

m c (0|Ci ( Si)|0) 

of 0((v 2 )) to be about 0.1 from Eq.( )25l) or Eq. f l27l) . which is about 2 ~ 3 times smaller than 
the result calculated within a potential model [3 li, l34j or within the QCD sum rules 35]. Note 
that, differing from theirs, our matrix elements are extracted by fitting experimental data, 
which depend on the chosen processes and experimental errors. Moreover, the short- distance 
coefficient of the v 2 correction term in Eq. fl23b ) for J/ip — >■ LH is quite large, implying that 
the corresponding long-distance matrix element could be rather small. These of course will 
have uncertainties, compared with other methods for calculating the long-distance matrix 
elements. Nevertheless, our method using experimental data to extract the matrix elements, 
provides an independent and self-consistent way to determine the matrix elements. So in 
this work we will use these experimentally extracted long-distance matrix elements to give 
numerical predictions. 
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C. Numerical Result 



The expression of Nq and N\ are too complicated to be shown here and we only give 
the analytical expression of the differential cross section. With J PC conservation and gauge 
invariance, the general form of the differential cross section of unpolarized J/ij) production 
in e + e~ annihilation can be expressed as 
d 2 cr 

— — — (e + e~ 7* -> i> + gg) = S ( Zl )(l + a («i)s?) + S v { Zl ){l + a^x?), (29) 

(XZ\(IX\ 

where the first term on the right-hand side of Eq. [29] is the LO result in v c , and the second 
one is the relativistic correction term. Our results of 5o(zi) and ao(zi) are in consistency 
with those in And the expressions of S v (zi) and S v (zi)ct v (zi) are too complicated and 
they will not be given here. 

Setting a s = 0.26, \fs = 10.6GeV and m c = 1.5GeV and integrating out X\ and Z\ 
numerically, we find when choosing the values of the matrix elements in Eq. [25], the LO 
result is 

<Tio(e + + e~ -»■ J/tp + gg) = 202 fb. (30) 

and the relativistic correction is 55 fb, which gives about 27% enhancement, then resulting 
in the NLO result in v 2 : 

<TNLO(t,g)(e + + e" ->• J /iff + gg) = 257 fb. (31) 
If we choose the values of the matrix elements in Eq. [27] . the NLO result in v 2 becomes 

<tnlo(„2) ( e+ + e_ + 99) = 480 fb, (32) 

and the relativistic correction enhancement is 1.22. When the charm-quark mass varies from 
1.4GeV to 1.6GeV, the LO and NLO cross sections as function of m c are shown in Fig. [2] 
with the long-distance matrix elements, respectively, in Eq. fl26j) and Eq. (f28|) . As mentioned 
above, we do not expand the three-body phase space by v 2 and assume Mj/^ = 2m c for 
simplicity. Instead, if we used the physical mass of Mj/^p = 3.097GeV, which includes v 2 
kinematic and binding energy corrections, when doing the phase space integrals, the results 
given in Eq. (131 j) and Eq. (l32l turned to be 253fb and 475fb, respectively. Comparing those 
results, one can see that the differences due to replacing Mj/^ = 3.097GeV with Mj/^ = 2m c 
in the phase space integration can be neglected. 
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FIG. 2: Relativistic corrections to c(e + e~ — > J/ip + gg) as functions of m c with long-distance 
matrix elements determined from Eq. (|26p (left-panel) and Eq. (j28j) (right-panel) . In each figure the 
lower curve is the LO result and the upper one is the result including v 2 c corrections. 



To be consistent, when we calculate the production cross section including both the 0(v 2 ) 
relativistic correction and the 0(a s ) radiative correction, we should first use the relativistic 
correction results obtained by adopting Eqs. ()27ll28l) . where the 0(a s ) radiative corrections 
to the matrix elements are included, and then further include the 0(a s ) radiative corrections 
to the production short-distance coefficients. This will be discussed in next section. 



IV. DISCUSSION AND SUMMARY 



In above the order v 2 relativistic effect is considered for the color-singlet e + e _ — > 
J/i/j + gg process in the framework of NRQCD. We calculate the short- distance coeffi- 
cients perturbatively and find the ratio of the coefficient of (0\Vf ( 3 <Si)|0) / m 2 c to that of 
(0\Of i 3 Si)\0) is about 2.5. Denoting the ratio of the long-distance matrix elements by 
(0\Pf( 3 Si)\0)/(m 2 c (0\Of( 3 Si)\0)) = <^> [ai], then the enhancement of the relativistic cor- 
rection can be expressed by 2.5(t>f). Using the matrix elements given in Eq.(25j and Eq. [27] 
as inputs, we predict the enhancements of relativistic corrections are 22% and 27%, respec- 
tively, which are as important as the NLO QCD corrections [18]. It can also be found that 
when including the relativistic corrections the m c dependence is also improved a little. If we 
determine the matrix element from the J/ip — >■ e + e~ process in Eq.|23h] without including 
the relativistic and QCD corrections, we get \ — — = 0.251GeV and the LO result of 



11 



e + e~ — > J/ip + gg is only about 173fb. Then we can see that the relativistic corrections can 
enhance both the short- distance coefficients and the long-distance matrix elements. 

Next, we further include the 0(a s ) radiative corrections to the production short-distance 
coefficients. In Ref.[18||. the authors obtain that with a s (2m c ) = 0.259 the K-factor of NLO 
QCD result to LO QCD result is 1.20. Then using the matrix elements in Eq. [27], which 
also include the NLO QCD corrections in decay processes, we find that after including the 
QCD corrections [l^l the combined cross section onlooi,^) ^ or mc = l-5GeV is 

o-NLO(^,a 3 )(e + e- ->■ J/V> + gg) ^ 480/1.22 x (1 + 0.22 + 0.20) ~ 560fb, (33) 

where on the right-hand side of Eq. (1331) the number 480 fb comes from Eq. (l32|) . and in the 
summation the enhancement factor 0.22 is due to relativistic correction, while the enhance- 
ment factor 0.20 due to QCD radiative correction. Note that all the above contributions 
come from the color-singlet part. The LO color-octet contribution of e + e~ — > J/ip + g can 



be estimated as large as 0.27pb 11], but this apparently depends on the chosen values of the 
color-octet matrix elements. If using this estimate [11] for the color-octet contribution, then 
the prediction of NRQCD for the J/ip + X non _ cE cross section at £>-factories would become 
0.83pb, which is almost twice as large as the measured central value 0.43pb by Belle [ijj]. 



In fact, from Eq. fl33|) we see that after including the QCD and relativistic corrections the 
color-singlet contribution alone has saturated the measured value of J/tj) + X non _ c5 cross 
section, and thus there seems no need for the color-octet contribution. However, we must 
pay attention to possible uncertainties before we can draw a firm conclusion. First, on the 
experimental side, there is a large uncertainty of Belle's result in the pj/^p < 2.0GeV region, 
and furthermore the total cross sections of inclusive J/ip production measured by BaBar, 
Belle, and CLEO are not consistent. Second, on the theoretical side, we should take into 
account the uncertainty due to the choice of renormalization scale /i in the calculation of the 
short- distance coefficients. In doing the latter, we choose the largest value fi = and 



m c = 1.5GeV, and find the K-factor of NLO QCD correction to be 1.48 (see Ref.[l8|). We 
then combine the relativistic correction with the QCD radiative correction, and find that at 
NLO in v\ and ot s the total cross section of direct J/0 production becomes 

o"NLO(i,2, as )(e + e" ->• J/ip + gg) ~ 437fb (34) 

for /i = y/s/2 and m c = 1.5GeV. We see that although the cross section is decreased as 
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compared with that for /x = 2m c and m c = 1.5GeV, the predicted cross section is still a 
little larger than the central value of the latest Belle result. 

Moreover, the cross sections obtained in Eq. (13"3~j) and Eq. flM]) are the direct J/ij) produc- 
tion rates, not including the feed-down contribution from higher charmonium states. If the 
feed-down contribution is included, the prompt J lib production cross section, which is the 

nn 

measured value by Belle, will be further enhanced by a factor of about 1.3 (see Refs.[18|, [19| 
for discussions on the feed-down contribution). Therefore, the theoretical cross section for 
the J fib prompt production calculated at NLO in a. and v 2 in NRQCD will exceed or sat- 
urate the latest value in Eq.flS]) observed by Belle[l4j, despite of theoretical uncertainties 
related to the choice of input parameters, e.g., /i, m c , and the color-singlet matrix elements. 

In summary, we find the 0(v 2 ) relativistic correction to enhance the cross section of 
J fib production in the color-singlet process e + e~ — >• J '/ip + gg by a factor of 20-30%, which 
is comparable to the enhancement caused by the 0(a s ) radiative correction 18l. [19| . As 
the consequence of including both the 0(a s ) radiative correction and the 0(v 2 ) relativistic 
correction, the color-singlet contribution to e + e~ — > J jib + gg has saturated the latest 
observed cross section by Belle for e + e _ — > J/ib + X non _ cE at 5-factories, thus leaving little 
room for the color-octet contribution. This gives a very stringent constraint on the color- 
octet contribution, and may imply that the values of color-octet matrix elements are much 
smaller than expected earlier by using the naive velocity scaling rules or extracted from fitting 
experimental data with the LO results. To reduce the theoretical uncertainties, further 
investigations for the higher order (both in a s and v 2 ) corrections are needed. Moreover, 
comparisons between various experimental measurements and theoretical predictions are 
certainly helpful to clarify this important issue concerning the color-octet mechanism. 
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